Renormalization group study of a kinetically constrained model for strong glasses 
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We derive a dynamic field theory for a kinetically constrained model, based on the Fredrickson- 
Andersen model, which we expect to describe the properties of an Arrhenius (strong) supercooled 
liquid at the coarse-grained level. We study this field theory using the renormalization group. For 
mesoscopic length and time scales, and for space dimension d > 2, the behaviour of the model 
is governed by a zero-temperature dynamical critical point in the directed percolation universality 
class. We argue that in d — 1 its behaviour is that of compact directed percolation. We perform 
detailed numerical simulations of the corresponding Fredrickson-Andersen model on the lattice in 
various dimensions, and find reasonable quantitative agreement with the field theory predictions. 



PACS numbers: 64.60.Cn, 47.20.Bp, 47.54.+r, 05.45.-a 



I. INTRODUCTION 

The discovery that supercooled liquids, whose struc- 
tures are essentially homogeneous and featureless, are dy- 
namically highly heterogeneous is arguably the most im- 
portant recent development in the long-standing problem 
of the glass transition 0, 0, Q • Dynamic heterogeneity 
has been observed both experimentally in deeply super- 
cooled liquids Q , and in numerical simulations of mildly 
supercooled liquids 0,0- In addition, dynamic hetero- 
geneity has been observed experimentally in colloidal sus- 
pensions • For recent reviews see |E IE El El • 

Understanding dynamic heterogeneity is a crucial step 
towards understanding the glass transition. Dynamic 
heterogeneity implies that the slow dynamics of glass- 
formers is dominated by spatial fluctuations, a feature 
discarded from the start in homogeneous approaches like 
mode coupling theories Hal or other mean-field treat- 
ments EI; see however |l4 ITU Il6l El- Moreover, 
the absence of growing structural length scales has been 
the main obstacle to the application to supercooled liq- 
uids of many of the tools used so successfully to ana- 
lyze conventional phase transitions. The existence of dy- 
namic heterogeneity, however, implies that the increase 
in timescales as the glass transition is approached is as- 
sociated with growing length scales of d y namica lly, not 
statically, correlated regions of space 

Ef (Refs. El H3 

offer alternative thermodynamic viewpoints). This sug- 
gests that supercooled liquids might display universal dy- 
namical scaling, by analogy with conventional dynamical 
critical phenomena. This scaling behaviour could then 
be studied by standard renormalization group (RG) tech- 
niques. This is the issue we address in detail in this paper, 
which is a follow-up to our recent Letter [2l| . 

Our starting point is the coarse-grained real-space de- 
scription of glass- formers developed in 0, I22L 123 ] which 
places dynamic heterogeneity at its core. It is a meso- 
scopic approach, based on two observations. First, at 
low temperature very few particles in a supercooled liq- 
uid are mobile, and these mobility excitations are local- 



ized in space; and second, mobile regions of the liquid 
are needed to allow neighbouring regions to themselves 
become mobile. This second observation is the concept 
of dynamic facilitation (23, HE HE H3- This picture of 
supercooled liquids can be straightforwardly cast as a 
dynamical field theory and its scaling behaviour deter- 
mined via RG. We find that, generically, scaling prop- 
erties are governed by a zero-temperature critical point. 
We study in detail the simplest case, that of isotropic dy- 
namic facilitation, which we expect to model an Arrhe- 
nius, or strong, glass-former. We show that its low tem- 
perature dynamics is controlled by a zero-temperature 
critical point, which for dimensions d > 2 is that of di- 
rected percolation (DP) |2g|. We argue that in d = 1 
the model belongs to the universality class of compact 
directed percolation (CDP) [2^. Our theoretical predic- 
tions compare favourably with our results from numerical 
simulations of the Fredrickson-Andersen (FA) model [25| , 
the lattice model on which the field theory is based. 

This paper is organized as follows. We derive a field 
theoretic description of a generic system possessing con- 
strained dynamics in Section II, and we discuss in Sec- 
tion III the physical interpretation of the field theory for 
the special case of an isotropically-constrained model. In 
Section IV we study this field theory for d > 2 using RG, 
and in Section V we discuss the special case of d = 1. In 
Section VI we compare the theoretical predictions to sim- 
ulations of the FA model in various dimensions. Section 
VII contains a summary of our results and conclusions. 



II. DERIVATION OF THE FIELD THEORY 

We build an effective model for glass-formers as fol- 
lows (2^- We coarse-grain a supercooled fluid in d spa- 
tial dimensions into cells of linear size of the order of the 
static correlation length, as given by the pair correlation 
function. We assign to each cell a scalar mobility, n*, 
whose value is chosen by further coarse-graining the sys- 
tem over a microscopic time scale. Mobile regions carry a 
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free energy cost, and when mobility is low we do not ex- 
pect interactions between cells to be important. Adopt- 
ing a thermal language, we expect static equilibrium to 
be determined by a non- interacting Hamiltonian |27|. 

N 

H = J2m- (1) 

i=l 

At low mobility, the distinction between single and mul- 
tiple occupancy is probably unimportant, and we assume 
the latter case for technical simplicity. The question of 
whether the field theories for versions of a system with 
single ('fermionic') or multiple ('bosonic') occupancies lie 
in the same universality class is unresolved 29] , but if the 
distinction matters it is likely to matter more in d = 1 
than in higher dimensions. We shall ignore this subtlety. 

We define the dynamics of the mobility field by a mas- 
ter equation, 

d t P{n,t) = Y J d{n) AP(n,i), (2) 

i 

where P (n, t) is the probability that the system has con- 
figuration n = {rii} at time t. Equation (01 shows clearly 
the two ingredients of our model. The first, the existence 
of local quanta of mobility, is encoded by the local op- 
erators Ci. For non-conserved dynamics we choose these 
to describe creation and destruction of mobility at site i, 

tiP{m,t) = j(ni + l)P(ni + l,t) + P P(ni-l,t) 
-(- f n i + p)P(n i ,t), (3) 

where the dependence of P on cells other than i has been 
suppressed. The rates for mobility destruction, 7, and 
creation, p, are chosen so that © obeys detailed bal- 
ance with respect to at low temperature. This means 
that the stationary solution of the master equation must 
equal the Gibbs distribution. Equation Q gives rise to 
the Gibbs distribution P cq (n) = ^(1 - e - 1 / T ) e -"*/ T , 
whereas the master equation (|2J) has the stationary solu- 
tion 

P(n,i^oo) = ne- p/T (^" , ^y. (4) 

Equation (0J will reduce to P C q(n) provided p/j = e~ x / T 
and Thus detailed balance with respect to 

holds only at low temperature, and we will hereafter as- 
sume that T « 1. For convenience we write p/j w c, 
where c = (rii); angle brackets (• • • ) denote an equilib- 
rium, or thermal, average. The thermal concentration of 
excitations c is the control parameter of the model. 

The second ingredient of our model is the kinetic con- 
straint, Ci {{n}), which must suppress the dynamics of 
cell i if surrounded by immobile regions. It cannot de- 
pend on rii itself if J5J) is to satisfy detailed balance. To 
reflect the local nature of dynamic facilitation we allow 
Ci to depend only on the nearest neighbours of i [23] and 
require that Ci is small when local mobility is scarce. 



One can derive the large time and length scale be- 
haviour of the model defined by Eqs. from an 
analysis of the corresponding field theory. The technique 
to recast the master equation © as a field theory is 
standard [3(], • O ne introduces a set of bosonic cre- 
ation and annihilation operators for each site i, a\ and 
<H, satisfying = and defines a set of states 
\n) = (a^™ |0), such that 

a] \m) = I ra< + 1) , a, \rii) = n % \m - 1} . (5) 

The vacuum ket |0) is defined by a |0) = 0. One passes 
to a Fock space via a state vector 

i*(i))^^p(n,t)n«r i -io). (6) 

{«.} i 

The master equation J3J then assumes the form of a Eu- 
clidean Schrodinger equation, 

$!*(*)>= -A|*(f)), (7) 

with H = ^2iCi({ojaj})H^ . The unconstrained piece 
reads 

H? ] = - 1 {a l -a\a l )-p{a\~l), (8) 

which describes the creation and destruction of bosonic 
excitations with rates p and 7. The evolution operator 
e~ Ht can then be represented as a coherent state path 
integral weighted by the dynamical action [3l| 

s[<t>*,<t>iM = Y<j° dt + W^)) , (9) 

where we have suppressed boundary terms coming from 
the system's initial state vector. The fields 4>*(t) and 
4>i{i) are the complex surrogates of a\ and a.;, respec- 
tively, but must now be treated as independent fields 
and not complex conjugates. The Hamiltonian Hi has 
the same functional form as JSJ with the bosonic opera- 
tors replaced by the complex fields. At the level of the 
first moment we have (n,) = (</>j), and so we may re- 
gard (f>i as a complex mobility field. Higher moments of 
rii and <pi are not so simply related, however: for ex- 
ample, (n|) = {4>f) + (fa) 32]. The last step in the 
passage to a field theory is to take the continuum limit, 
according to J2i ~~* a <l I d d x, <Pi(t) — ► a d 0(x, t), and 
(j)*(t) —¥ </>*(x, t), where a is the lattice parameter. 

The definition of the model is completed by specifying 
the functional form of the kinetic constraint. The sim- 
plest non-trivial form is an isotropic facilitation function, 
Ci = rij , where the sum is over nearest neighbours 
of site i. With this choice we expect our model to be 
in the same universality class as the one-spin facilitated 
Fredrickson- Andersen model in d dimensions [2^, |27| . 
Different choices for the operators and C lead to 

field theoretical versions of more complicated facilitated 



models. A diffusive H^°\ for example, would correspond 
to a constrained lattice gas like that of Kob and Ander- 
sen [13, El; an asymmetric C to the East model p7ll83| 
and its generalizations p^ . 

In the continuum limit the isotropic constraint reads 

£>*0i « (2d + a 2 V 2 + ---)0*^, (10) 

3 

where ' ■ • • ' denotes higher-order gradient terms irrelevant 
in the RG sense in the long time and wavelength limit. 
Terms linear in the spatial gradient vanish because the 
constraint is isotropic. Consequently, the dynamics of 
the model is nearly diffusive, perturbed by fluctuations 
in low dimensions. 

To derive the dynamic action it is convenient to make 
a linear shift of the response field, tp* — > 1 + [3l| . 
in the Hamiltonian. This is done for the following rea- 
son. Expectation values in this formalism are given by 
by (A) = {s\A\Mt\u where = e~ At |tf(0)) and 

(s| = (0| e^i ai |3l|,|32|. The projection state (s| is intro- 
duced because the usual quantum mechanical expression 
\A\ is bilinear in the probability P. If one 

wishes to apply Wick's theorem, one must commute the 
factor e^i ai to the right hand side of the bracket; the 
consequent shift <fi* — > 1 + <p in the Hamiltonian follows 
from the identity e a f(a\ a) = /(l + a\ a)e a , and corre- 
sponds to a change of integration variables. It therefore 
does not change the properties of the system under renor- 
malization. However, it can obscure important symme- 
tries of the model in question, and so should be made 
with care |3l| . 

The dynamic action now follows from Equations JH}, 
@ and i|l(J|) . suitably shifted, and reads 

S$, 0, t ] = j d d x J °dt[$(dt- D V 2 - r ) <j> 

+ i/W V 2 )0 + 0>(A O 2) + z/£ 2) V 2 )# 
-^o + ctoV 2 )^) . (11) 

We have defined Ag 1 ' 2 * 1 = 2da d j, r$ = Vo = 2dp, 
j/q 1 ' 2 ^ = 7a d+2 and cr = a 2 p. We write D = er to 
emphasise the emergence of a diffusive term, although in 
the unshifted model there is no purely diffusive process. 
We have omitted higher-order gradient terms, and sup- 
pressed boundary contributions coming from initial and 
projection states. Equation (fTT|) is the starting point for 
our RG analysis. 
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FIG. 1: Elements of the simplified action 1121 . From left 
to right: branching and coagulation vertices, and diffusive 
propagator. Time runs from left to right. These processes 
appear in the space-time trajectories of simulations of the 
lattice FA model in Fig. |5] 



of the model by dropping all but the most relevant terms 
from the action to give (see Section llV|) 

S = J d d xdt 4> (d t - D V 2 - r ) <t> 

+ J d d xdt (A o 00 2 - v 4) 2 (j)) . (12) 

The first term is the bare propagator of the theory, the 
renormalized version of which corresponds to the prob- 
ability that two sites separated in space and time are 
connected by an unbroken chain of mobile sites ■ The 
second and third terms are the vertices corresponding to 
coagulation and branching interactions, respectively. In 
the usual way [jjifj one associates with each term in the 
action a diagram, as in Fig. ^ 

The physical processes corresponding to the terms in 
(|12|l or the diagrams in Fig. ^ can be seen in numerical 
simulations. In Fig. [2] we show a typical space-time tra- 
jectory 0| for the FA model in 1 + 1 dimensions. The 
wandering of excitation lines corresponds to the diffusion 
of isolated defects. Diffusion appears in the propagator 
as a result of the shift (p* — > 1 + cf> applied to the term 
</>*V 2 0*</>, which enters Eq. 0. This term corresponds to 
nearest-neighbour facilitated mobility creation with rate 
oc c, and so diffusion in our model results from facilitated 
creation (branching), followed by facilitated destruction 

(coagulation): {T ^TT T} ~ c0V 2 0. 

Branching and coagulation events can be clearly iden- 
tified: one of the latter is enlarged in the lower left of the 
figure. These events correspond to fluctuations. In low 
dimensions, where fluctuations are important, branching 
and coagulation events renormalize the bare propagator 
of the theory, meaning that excitation lines joining two 
sites are dressed by bubbles. In low dimensions one must 
therefore resort to RG in order to account for fluctuation 
effects in a controlled way. 



III. PHYSICAL INTERPRETATION OF THE 
ACTION 

Equation Ijlll) has the form of an action for a sin- 
gle species branching and coagulation diffusion-limited 
reaction psl l35l | with additional momentum-dependent 
terms. We can see how this action governs the behaviour 



IV. RG ANALYSIS OF THE ACTION 

A. Langevin equation of motion 

By making stationary variations of the action (| 1 1 1> with 
respect to the response field 4>, 5S/8<j> — 0, we obtain the 




FIG. 2: A space-time trajectory for the FA model in 1 + 1 
dimensions; time runs horizontally from left to right; space 
vertically. Mobile sites are black. The events corresponding 
to the diagrams in Fig.^can be clearly seen. The wandering 
of the excitation lines corresponds to diffusion; branching and 
coagulation events act to renormalize the excitation lines. In 
the lower left a coagulation event is shown enlarged. 



Langevin equation of motion for the field <f>: 

d t <t>{t) = DoV 2 ^ + r o - A o 2 + tj{x, t), (13) 
where the noise rj satisfies (r)(x,t)) = 0, and 

(r)(x, t)r,(x',t')) = 2 (v <f> - A< V) 8(x - x')S(t - t'). 

(14) 

We neglect diffusive noise. The noise-noise correlator 
(|14fl describes stochastic fluctuations of the mobility field 
(j), and comes from the coefficient of the terms in the ac- 
tion quadratic in (j>. We see that <fb4"|) describes a competi- 
tion between mobility correlations and anti-correlations, 
induced by branching and coagulation, respectively. If 
for example a branching event occurs, a particle will find 
itself with more nearest neighbours than one would ex- 
pect from a mean-field argument. In the long time and 
wavelength limit, for d > 2, we show below that the sec- 
ond term in Eq. 114|) is irrelevant, and may be dropped. 
Equations (|13[1 and 114|l then constitute the well-known 
Langevin equation for DP [28j , albeit with a positive def- 
inite mass term. 

The mean-field approximation consists of dropping the 
noise and diffusion terms from (|13[) . The resulting equa- 
tion possesses a dynamic critical point at c oc tq = 0, 
or T = 0. For T > 0, in the non-equilibrium regime, 
the density cf> approaches its thermal expectation value 
c = p/7 exponentially quickly. At T — the decay 
becomes algebraic, (f>(t — » 00) ~ (\ot)~ . Whether at 
equilibrium or not, the mean-field equation admits the 
critical exponents i>\\ = 1 and (3=1. The former de- 
scribes the growth of time scales £|| near criticality, via 
£|| ~ c _!y n ; the latter is the order parameter exponent, 
defined as the long-time scaling of the density in terms 
of the control parameter, n(t — > 00) ~ c 13 . By restoring 
the diffusive term, the spatial exponent v\_ = 
analogously to ^|| , may be identified. 

In the following section we show that fluctuations al- 
ter these predictions in low dimensions, by virtue of 



^, defined 



(a) 



(b) 




(c) 



FIG. 3: Vertices corresponding to terms in the action I1H . 
(a) From left to right, diagrams corresponding to the vertices 



-A 



-A 



(2) 



and i>q gi. 



The dot denotes a momentum 



dependence q\ on one of the incoming legs, (b) The structure 
of the lowest-order corrections to the propagator, showing how 
the effective couplings emerge. The coupling shown here is x. 
Similar diagrams with dots on right, left and both vertices 
correspond to propagator renormalization controlled by the 
couplings y, z and u, respectively, (c) Diagrams renormal- 
izing the coupling x. Because these diagrams are symmetric 
with respect to incoming and outgoing momenta and frequen- 
cies, the contribution from each is doubled. Note that (b) also 
contributes to the renormalization of x, via the renormaliza- 
tion of D. 



endowing space and time scale exponents with small 
dimension-dependent corrections. The exponent /3, how- 
ever, remains unchanged. We argue that because our 
model possesses detailed balance, which ensures that 
(4>(t — > 00)) = (4>) = c, this fixes (3 to unity. This may 
also be inferred from the invariance of the unshifted ac- 
tion under the transformation (6^6*) — > (c</>* , c -1 (/>) , and 
the consequent Ward identity j37| ■ 



B. Dimensional analysis 

We identify the upper critical dimension of the model 
via dimensional analysis [3(| |3^| . We rescale space ac- 
cording to x — > xp 1 / 2 , in order to remove the temperature 
dependence from the diffusion coefficient. Note that this 
rescaling is not valid at c = 0. The action Hll|) then reads 

as before, with rescaled parameters Aq 1,2 -* = 2da d p 



7, 

r = v = 2dp, <7q = a 2 , and J^ 1 ' 2 ^ = 7/J~ 1 ~^a d+2 . We 
show some of the diagrams corresponding to these cou- 
plings in Fig. 

To perform a scaling analysis, we identify the effective 
couplings emerging from the action. These follow from 
the structure of the diagram shown in Fig. Et>, and are 

_ V A _ cr A _ V V Q _ (J h> Q 

X V — — Tv) ) Vo = — T^T~ ' z = — 7-^9 ' u 



n 2 



n 2 



D 2 



n 2 



(15) 

The factors of Do come from the explicit evaluation of 

the integrals associated with the diagrams. To H15|) we 

(2) (2) 
add go = A /Do and ho = v /Do, which couple to 

four-point vertices: see Fig. [3] Dimensional analysis re- 
veals that the upper critical dimension is 4, at which the 
most relevant coupling, xo, is marginal. Renormalization 
effects must therefore be taken into account for d < 4. 
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Above d — 4 the classical (fluctuation-free) predictions 
apply. Other couplings become relevant below d = 2, 
and we shall therefore restrict our analysis to 2 < d < 4. 
Dimension d = 1 is treated separately in Section Ivl 



C. DP fixed point 

We employ the usual field-theoretic renormalization 
group scheme using dimensional regularization 

in d = 4 — e dimensions to identify the unphysical ultra- 
violet (UV: short time and distance) poles of the ver- 
tex functions r^'^ of the theory. The vertex func- 
tions r(^ w ) consist of all one-particle-irreducible dia- 
grams with N outgoing and N incoming amputated lines. 
Their UV poles result from exchanging a lattice model, 
which is regularised at short distances, for a continuum 
field theory, which is not. But by invoking universality, 
which says that the behaviour of a system approaching 
criticality is governed by a small number of relevant pa- 
rameters, we recognise that the UV poles correspond to 
irrelevant microscopic degrees of freedom. By removing 
these poles we both render our theory finite, and, via 
scale-invariance and dimensional analysis, infer its phys- 
ically important infra-red (IR: large time and distance) 
scaling |38j. We shall work to one-loop order, and use 
dimensional regularization and minimal subtraction |36| . 

We introduce the following renormalized counterparts 
of the fields and couplings appearing in 



V = Z^ZyVQfl 



R ' 

-2 

) 

D 



r = Z^,Z r (r - r 0c )n 2 , 
Z^ZdDq, 



where r$ c is the additive counterterm introduced to cure 
the quadratic divergence of the vertex function I^ 1 ' 1 ). 
We have introduced an arbitrary momentum scale /i in 
order to render the couplings dimensionless, and have 
chosen to allocate dimensions to the fields according to 
[0] = 1, [4>] = The predictions of the theory must be 
independent of this allocation. 

We define the multiplicative renormalization factors Z 
as follows. From the propagator couplings, we fix mass, 
field and diffusion constant renormalization via 



NP 



= 1, 



NP 



-r^ 1) (0,0) = r/x 



D. 

2 



(16) 



while for the couplings comprising xq we impose the con- 
ditions 



A« = i it 1 - 2 ) 



,6,-1 



NP 



NP 



(17) 



The subscript NP stands for 'normalization point', and 
is the value of the external momentum scale at which 



we evaluate the vertex functions. It can be chosen for 
convenience, provided that it lies outside the IR-singular 
region; we take NP = (iuj,q 2 ,r) = (2D^ 2 ,0,0). Note 
that this choice corresponds to the system at criticality, 
which for finite T is an approximation. For non-zero T 
one must retain the mass term in the propagator. This 
leads to the emergence of an effective coupling that flows 
logarithmically to zero, signaling a crossover to a mas- 
sive, classical fixed point. We will discuss this case in 
Section ITVFl 

We first assume that for 2 < d < 4 the couplings other 
than xq are irrelevant, and hence the action l|llfl reduces 
to that of DP (we shall call x the 'DP coupling'). We shall 
find that those couplings which are marginal in d = 2 at 
the classical fixed point are rendered irrelevant at the 
DP fixed point. Hence we expect to see DP scaling for 
2 < d < 4. We find, to one-loop order, the well-known 
Z-factors 0, 

v 1 , Ad - 
Z r = l- "^-^oA 1 , 



where A d = 4(47r)- d / 2 r(3 - d/2). Note that the cu- 
bic vertices renormalize identically as a consequence of 
a Ward identity. The renormalization factor associated 
with x = X^v/D 2 is therefore Z x = Z^,Z xW Z V Z D 2 
1 — (3^4ci/2e)a;o/i _c + 0(e 2 ). Insofar as one can ignore 
the propagator mass, the rescaled coupling x — > A^x = 
AdZxXQpT^ changes with the observation scale /i accord- 
ing to 



Z D 


= 1 - 


A d _ 
8e 




= 1 - 


A d _ 
— xqV- 
e 


Z v 


= 1 - 


A d _ 

— xoV 
e 



(3 X 



/' 



dx 



(19) 



If we parameterize the change in the observation scale by 
fj, —> = fit, we can solve (|T§)) for x(l): 



x(£) 



t + 1 



(20) 



Thus x — > x* = 2e/3 as £ — > 0, because e > 0. Since I w 1 
and < < 1 correspond respectively to microscopic and 
macroscopic length and time scales, x* is an IR-stable 
fixed point. At this fixed point the critical exponents of 
the theory are independent of its microscopic parameters, 
and so are 'universal'. We therefore expect the model to 
display scaling behaviour independent of its microscopic 
details for very low temperatures. This scaling behaviour 
belongs to the universality class of directed percolation. 

Having assumed the non-DP couplings in the action 
are irrelevant for 2 < d < 4, we shall now justify this 
assumption. These couplings are indeed irrelevant at the 
classical fixed point, as one can verify from l|llfl by di- 
mensional analysis. We find that they remain irrelevant 
at the DP fixed point. Further, those couplings (y,z,g) 
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which are marginal in d = 2 at the classical fixed point 
are rendered irrelevant at the DP fixed point. Hence we 
expect to see DP scaling in d = 2, also. Defining a renor- 
malization scheme in a similar manner to before, 



_ 7 ,,a 



,(2) _ 



= Z?Z„ (a)I /<V, (21) 



where 



A< 2 



-r 



(2,2) 



OX 



NP 

a=d q2 T%' 1 \q,u) 

-(1,2), 



v^=d q *T^>{q,u) 



NP 



NP 



NP 



(22) 



we find to one-loop order 



7 - 7 - 1 2Ad - 

e 

Z A (2) =1 X /i 

e 



(23) 



The corrections to the Gaussian scaling dimensions of 
these couplings may then be calculated. The correc- 
tion to the Gaussian eigenvalue of g = X^'/D is de- 



termined by Z„ 



Z^Z 



1 - 5ac /i~7(8e). So 



7* = lidpki(g/go)\ x=xi , = d - 2 + 5e/12. Thus 5 is less 
relevant at the DP fixed point than at the Gaussian fixed 
point, and for d > 2 may safely be ignored. So too may 
h. In a similar way we find that 7* = 7* = d — 2 + 5e/3, 
and 7* = d + 7e/3, all of which are irrelevant for d > 2 
at the DP fixed point. Thus for 2 < d < 4 the scaling 
properties of the model near criticality are those of DP. 

The critical exponents of the model then follow from 
standard arguments [3 |3(|. They are, to 0(e), 

(v®* 3 , = (| + 1 + j^)- The temporal exponent 

appropriate for comparing these predictions with numer- 
ical simulations is 1 + i/u . The additional factor of unity 
arises because the microscopic timescale associated with 
the system goes itself as ~ c. The non-trivial value of j3, 
f3 DP = 1 — e/6, cannot be observed for a model such as 
ours which eventually equilibrates, as discussed above. 



D. Critical temperature 

In general, systems in the DP universality class, such 
as directed bond percolation in 1 + 1 dimensions, exhibit 
a continuous phase transition from an active to an ab- 
sorbing state at some finite value p c ^ of their control 
parameter p. Our model, for which p = c, displays no 
such transition. One can justify this difference on phys- 
ical grounds, as follows. If we interpret the mobility 
as the concentration of a chemical reactant A, then the 



KCM we study for d > 2 corresponds to a chemical re- 
action involving diffusion (A + *-* + A), branching 
(A + -> A + A) and coagulation (A + A -> A + 0). Re- 
call that the diffusive process arises from the mechanism 
of mobility creation facilitated by a nearest-neighbour 
site, and is made manifest only following a shift of the 
response field. DP corresponds to these three processes 
plus self destruction, A + ^> + 0. It is self destruc- 
tion that permits other systems in the DP universality 
class to undergo a phase transition at a finite value of 
the control parameter. Self destruction gives rise to a 
second-quantized operator 



H sd = -<Jo(ai - ojoi), 



(24) 



which, following a shift of the response field, results in a 
term in the action of the form ao4>4>. Thus the mean-field 
critical point becomes p c = uq. Near criticality, p c is in- 
creased above its mean-field value by fluctuations. This 
occurs because the DP noise-noise correlator is positive, 
and so coagulation is enhanced by the branching process: 
each particle finds itself with more neighbours with which 
it may coagulate than one would expect from a mean- 
field approximation. This enhanced coagulation enters 
the term which renormalizes the mass, effectively enhanc- 
ing self-destruction relative to branching, and shifting the 
critical percolation threshold upwards. 

Now self destruction is excluded by any dynamical rule 
preventing mobility destruction unless facilitated by a 
nearest neighbour. Moreover, no such process can be gen- 
erated under renormalization from only branching and 
coagulation processes whose respective rates are fixed by 
detailed balance. Hence we expect one-spin facilitated 
models in general to have a critical point at zero temper- 
ature. 

This argument may be made explicit for the model 
we study. By imposing the condition for criticality, — 
r^ 1,1 '^ = 0, q = 0,ro = ro c ), we find that, to one-loop 
order 



d-2 



roc 



(Xqvo/Dq) N d A 
l-(XoVo/D*)N d A d - 



(25) 



Here, for convenience, we have imposed an explicit 
wavevector cutoff A; the additive correction to the mass is 
formally equal to zero in dimensional regularization, and 
yet the physical shift of the critical temperature must 
be independent of the regularization scheme used |38| . 
We have introduced Nd{2ir)~ d {d — 4)~ 1 S d , where S d = 
27r d / 2 /r(<i/2) is the surface area of a d-dimensional hy- 
persphere, and N d = (d — A)N d / (d — 2). We also use the 
unsealed variables of the action (|llf> . in which Do oc c. 

From (|25|l we see that the critical bare mass changes 
sign as T — > from above. Ostensibly the critical tem- 
perature is then negative; physically, of course, it is zero. 
This is a consequence of the vanishing of fluctuations 
in the limit of zero temperature, which may be inferred 
from the vanishing in that limit of the branching vertex 
in the action. The diffusion term arises from nearest- 
neighbour-facilitated branching, and so must also vanish 
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in this limit. Thus there is no fluctuation-induced shift 
of the critical temperature which remains T c = 0. 

This is as we expect, if the field theory is a faithful rep- 
resentation of the original master equation. The master 
equation satisfies detailed balance at all temperatures, 
which means that it cannot admit an absorbing state: 
an absorbing state breaks detailed balance because it is 
a state that may be entered, but not left. Nonetheless, it 
is necessary to verify, as in Eq. (|25[1 , that there exists no 
finite-temperature absorbing state under coarse-graining 
of the master equation. The FA model, upon which the 
field theory is based, is known to have a critical point at 
zero temperature |27j ] . 



E. Crossover to classical behaviour 

For any T > the mass parameter tq oc c will be non- 
zero. Under renormalization, as discussed above, it will 
eventually become large, rendering our approximation of 
criticality incorrect. The system will thus for very large 
time and length scales exhibit classical scaling properties, 
with the associated simple exponents. 

We can quantify the emergence of the classical theory 
by retaining the mass term in the propagator 31]. If we 
write s = r/D, we find that 



dlnx{£) 

dint 
dins (I) 

d\nl 



(26) 
(27) 



where g(l) = x(£)(l + s(l)) d > 2 - 2 is an effective coupling. 
For small s, x would in the IR limit approach to the 
directed percolation fixed point. But s does not remain 
small, flowing as s(£) ~ £- 2 +o(0. if 

we introduce the 

scaled mass a = (l + s(l)~ 1 £ 2+0 ( c ^ , we find that in 
the large mass limit a — > 1 we obtain a logarithmically 
diminishing coupling, 



?M ~<?(i) 



1 + iV 4 



d)g(l)\n£- 



(28) 



The vanishing of the effective coupling signals the re- 
emergence of a classical theory: because g couples to 
diagrams renormalizing the propagator, its logarithmic 
vanishing results in a logarithmic crossover to classical 
exponents. 

Thus we should see DP scaling provided that temper- 
atures are small enough and time and length scales are 
not too large. The crossover temperature will be sys- 
tem dependent, because the prefactors of the flowing cou- 
plings are non-universal. For larger temperature or large 
enough length and time scales we expect to see a logarith- 
mic crossover to a classical theory. This is characterised, 
in the nonequilibrium regime, by exponential decay to 
the steady state, and in general by classical scaling be- 
haviour. 




time 



time 



FIG. 4: FA model in d = 1. The left panel shows an equilib- 
rium trajectory for the mobility field rii(t) atT — 0.3 (window 
size is L — 250 by At — 5000). The right panel shows the 
corresponding trajectory of the persistence field Pj(t); black 
sites denote those which have been or are mobile, and so sat- 
isfy Pi(t) = 0. The clusters generated by the dynamics of 
the Pi(t) are compact, and their scaling properties are that 
of CDP. 



V. DIMENSION d = 1 AND CDP 

For d = 1 DP scaling no longer holds. This is sig- 
naled in the field theory by the relevance of some of the 
non-DP couplings between d = 2 and d = 1, and the re- 
sulting profusion of uncontrollable singularities [3q|. In 
this section we argue that in d = 1 systems with single- 
spin isotropic facilitation, such as the FA model, belong 
instead to the universality class of compact directed per- 
colation (CDP) [28|. 

Consider the FA model in d = 1. The elementary 
order parameter of this model is the mobility field rii(t). 
Figure ^ (left panel) shows a portion of an equilibrium 
trajectory at T = 0.3. The connection between the FA 
model and CDP is made apparent by considering instead 
the corresponding persistence field Pi(t), i.e. the field 
which takes value Pi (t) = if site i has flipped by time 
t, and Pi(t) = 1 otherwise. The corresponding trajectory 
of Pi (t) in our example is shown in Fig. 0] (right panel) . 

Clearly, while the dynamics of rii(t) is reversible, that 
of Pi{t) is not. A related observation is that the clus- 
ters generated by the evolution of Pi(t) are compact, as 
seen in Fig. 21 The control parameter is again c, with 
c = corresponding to the transition between an active 
phase in which Pi{t) eventually becomes unity through- 
out the whole system, and an inactive phase, in which 
it does not. As before, the exponents and v±_ deter- 
mine the scaling of times, (Dt) ~ c _I/ H (with D k, c), 
and lengths, £ ~ cT v ^. Two further exponents deter- 
mine the asymptotic values of (Pi(t)). To extract these 
exponents it is convenient to define the transience func- 
tion Ti(t) = 1 — Pi(t): starting from an initial finite seed, 
(Ti(t — > oo)} ~ cr, with the dynamics running in the 
forward time direction; starting from a completely full 



lattice, (Ti(t 



-oo)) ~ c' , with the dynamics run- 
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ning backwards in time. Note that 0^/3' due to the 
irreversibility of Tj(i) (or P t {t)). 

The domains of Ti(t) spread only through diffusion of 
mobility excitations and interactions play no role. In 
this sense, the scaling behaviour of Ti{t) should be that 
of freely-diffusing domain walls, and coalescing domains. 
Examples of systems which behave similarly are the zero 
temperature Ising chain under Glauber dynamics [28j , or 
the reaction-diffusion system A + A — > Hi| . These 
indeed belong to the CDP universality class. 

CDP has the following exponents [28| : 



1 



2, vS 



1, 



CDP 



0, p 



/CDP 



1. (29) 



_FA _ , .FA 



These are precisely the values of the exponents of the 
FA model in d = 1. The time and length exponents 
= 1, giving the dynamic exponent 

/i/f* = 2 [H,|23|. Each site of a lattice which 
initially contains at least one excitation will eventually 
flip, and thus Ti(t — > oo) = 1 for all i, independently of 
c. We therefore have FA = 0, which is a consequence of 
ergodicity in the active phase. Conversely, if one takes 
a final state with all Tj = 1, and runs time backwards, 
the state at t — > — oo will have a density of excitations, 
and therefore of Ti , equal to c. This is a consequence of 
detailed balance. Hence 0' FA = 1 = ,CDP . 

We propose a field-theoretic justification for this be- 
haviour as follows. The Langevin equation of motion for 
4> is given by (|13[) and (|14H . At and above d — 2 the term 
in A 1 - 2 -* is irrelevant at the DP fixed point and may be 
dropped, leving us with the DP Langevin equation |28| . 
In d = 1, however, at the DP fixed point (assuming it 
exists), we have from our previous results the anoma- 
lous dimensions of the couplings appearing in the noise 
correlator: 



7a(2) =d-2 + e/3, 7* = 0. 



(30) 



We have calculated these dimensions using the prescrip- 
tion [(f)] = [(f)] = /i d / 2 , appropriate when the cubic vertices 
are considered independently. We see that A*- 2 -* and v are 
both marginal in d = 1. This is, we stress, a crude ap- 
proximation, because the calculation of the anomalous 
dimensions assumes the irrelevance of A^ 2 '. But it does 
suggest that here this assumption is inconsistent. Assum- 
ing that we can trust these exponents, there should then 
exist a fixed point controlled by \( 2 \ at which v is irrel- 
evant. Assuming this is so, and assuming further that 
the system can access this fixed point, this would leave 
the only vertices in the effective theory <jxf) 2 and (4><fi) , 
which allow no propagator renormalization. Hence z = 2 

exactly. For this theory the beta function is calculable 

(2) 

to all orders, since perturbation theory in X gives us 
a geometric series [39J. We then have a new fixed point, 
at which there exists a renormalized value of A^ 1 ' corre- 
sponding to an infinite value of its bare counterpart, Aq 1 ^ 



Thus A 



3> tqc/), giving the effective theory 




FIG. 5: Persistence function P(t) for the d = 3 FA model. 
From left to right: T = 5.0, 1.5, 1.0, 0.6, 0.4, 0.25, 0.17, 0.13, 
0.106 and 0.009. 



This is the Langevin equation for the CDP universality 
class [3^. We stress that this argument is conjecture 
only. A more rigorous analysis of the field theory would 
be required in order to justify this claim. 



VI. SIMULATIONS OF THE d = 3 FA MODEL 



The one-spin facilitated FA model |25j, |27| is the lat- 
tice model upon which the field theory of the previous 
sections is based. In this section we report the results of 
our large-scale numerical simulations of the equilibrium 
dynamics of the FA model in dimension d = 3, and com- 
pare these results to the predictions of the field theory. 
While the one dimensional FA model has been exten- 
sively studied by numerical simulations |27| . we are not 
aware of any detailed numerical study for d > 1. 

We consider the FA model on a cubic lattice with pe- 
riodic boundary conditions. The model is defined by the 
Hamiltonian and the isotropic dynamical rule 



Ci c 



Ci (l-c) 



(32) 



dtW) = A>V 2 - Af, 1 



'r)(x,t). 



(31) 



The kinetic contraint is Ci = 1 — Y[/j i) (1~ n j)i where (j, i) 
denotes nearest-neighbour pairs. We perform Monte- 
Carlo simulations of this model for several temperatures 
in the ran ge T £ [0.09, 5.0]. We use the continuous time 
algorithm |4(|, which is well-suited to this problem. The 
dynamical slow-down in this model is accompanied by 
the growth of a dynamic correlation length, and hence 
we must account for possible finite size effects. For in- 
stance, at T — 0.09 it was necessary (and perhaps even 
then not sufficient: see below) to use system sizes as large 
&sN = 160 3 . 
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FIG. 6: Arrhenius plot of the mean relaxation time in the 
d = 3 FA model. The DP exponent fits the data, while the 
classical exponent does not. 



A. Global dynamics 



< 2.5 




FIG. 7: Dimensionality dependence of the time exponent 
A = 1 + of the FA model, from numerical simulations in 
dimensions d = 1 to d = 6 (filled squares). The full line is the 
e expansion prediction. The values of the exponents agree 
within errorbars with those for DP for all d > 1, and and CDP 
for d — 1. The upper critical dimension of the FA model is 
£ A = d? P = 4. 



We first consider the spatially averaged dynamics. 
This may be probed via the mean persistence function, 




where Pi (t) is the single-site persistence function at time 
t, which takes value 1 if site i has not flipped up to time 
t, and value otherwise. Fig. [5]shows, as expected, that 
the dynamics slows down markedly when temperature is 
decreased below T k, 1.0, which marks the onset of slow 
dynamics in this model |23.l41|. 

We extract the mean relaxation time, t(T), via the 
usual relation P(t) = e _1 . The temperature dependence 
of t is shown in Fig. [ft] where various fits are also in- 
cluded. The high temperature behaviour is well described 
by a naive mean- field approximation 23], 

tmf-c- 1 . (34) 

This behaviour breaks down below T Q , where fluctuation- 
dominated dynamics becomes important. From our field 
theoretic arguments we expect that in the non-trivial 
scaling regime 

r~c~ A , A = 1 + w 2.1, (35) 

where the numerical value is the DP estimate in three 
dimensions [28j|. Fitting our data with the form r ~ c~ A 
we find 

A = 2.095 ±0.01, (36) 

as shown in Fig. El We include for comparison a fit using 
the Gaussian value of the exponent, A = 2, which is 
inconsistent with our data. 

We show in Fig. {7\ the results for similar simulations 
of the FA model in dimensions d from 1 to 6, together 



with the relevant DP exponent to 0(c 2 ), as tabulated in 
Ref. [2£|. We note that these results are consistent also 
with numerical simulations of systems in the DP univer- 
sality class j2^|. Our numerics also show that one-spin 
facilitated FA models display, in all dimensions, Arrhe- 
nius behaviour. They are thus coarse-grained models for 
strong glass- formers, as expected |22(. 

In summary, Fig. [7| strongly supports the RG predic- 
tion that the FA model exhibits non-classical scaling in 
low dimensions, consistent with DP behaviour for d > 2, 
and CDP behaviour in d = 1. 



B. Distribution of relaxation times 

The mean relaxation time t(T) captures only in part 
the relaxation behaviour of the model. We consider in 
this subsection the distribution of relaxation times, n(t), 
related to the mean persistence function via p3l | 

/oc 
dt'ir(t'). (37) 

These distributions are shown in Fig. |SJ 

A careful study of the functions Pit) and n(t) re- 
veals the following structure. At very large times, the 
persistence decays to in a purely exponential man- 
ner, P{t 3> t) ~ exp(— t/r). This is not the case in 
d = 1, where asymptotically the decay is described by a 
stretched exponential with stretching exponent (3 = 1/2. 
That stretched exponential behaviour is not seen in d = 3 
is consistent with the fact that strong glass-formers dis- 
play an almost-exponential relaxation pattern 0. 

Using t(T) as a unique fitting parameter does not al- 
low a satisfactory description of the whole decay of the 
persistence function: see Fig. |5| This figure shows that 
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FIG. 8: Distribution of relaxation times (symbols) and fits 
to Eq. JHHl (full lines) for temperatures T = 1.0, 0.4, 0.25, 
0.17, 0.13, 0.12, 0.106, and 0.09 (from left to right). 

there exists an 'additional short-time process', in the lan- 
guage of glass transition dynamical studies. 

Indeed, we find that fitting our data with the expres- 
sion 

?r(t) ~ t~ Q exp (-~) > ( 38 ) 

where a and r are free parameters, describes the distri- 
butions reasonably well over several decades: see Fig. [S] 
Often, data in the supercooled liquid literature are pre- 
sented in the frequency domain, because many decades 
can be accessed via e.g. dielectric spectroscopy j^. Fol- 
lowing this convention, we present frequency data ob- 
tained from the distribution of time scales via 

P(cj)= [ 7r(log(r))— — dlogr. (39) 

Interestingly, the short-time power law behaviour n(t) ~ 
t~ a observed in the distributions of time scales is also 
apparent in the frequency space as an 'additional process' 
on the high-frequency flank of the a relaxation, P"(cj) ~ 
(j-l+a. see Fig. |21 In this figure, the full lines correspond 
to fits of the main peak with a simple exponential, as 
discussed above. 

This feature is reminiscent of the 'high-frequency wing' 
discussed at length in the dielectric spectroscopy litera- 
ture 0- The wing is usually oberved in fragile glass- 
formers; unfortunately no dielectric data is available for 
strong glass-formers 43] . Other techniques, such as Pho- 
ton Correlation Spectroscopy, hint at the presence of an 
additional process in strong glass-formers similar to that 
observed in Fig. 44] . More experimental studies of 
the dynamics of strong glass-formers would be needed to 
confirm and quantify this similarity. 

C. Dynamic heterogeneity 

The growth of timescales in the FA model, r ~ c _A , is 
accompanied by growing spatial correlations, £ ~ c~ Ul - , 



t 




U) 



FIG. 9: Top: The fit of the persistence function with a simple 
exponential reveals an additional short-time process. Bottom: 
This is also true in Fourier space, where the additional process 
looks like a 'high-frequency wing'. 



as the system approaches its critical point at T = 0. 
These correlations are purely dynamical in origin, and 
give rise to dynamic heterogeneity Figure El il- 

lustrates this phenomenon in the FA model. We quantify 
the local dynamics via the persistence function Pi(t). For 
a given temperature we run the dynamics for a time t*, 
such that P(t*) = 1/2, meaning that half of the sites have 
flipped at least once. We colour white persistent (immo- 
bile) spins, for which Pi(t*) = 1, and black transient (cur- 
rently or previously mobile) spins, for which P%(t*) = 0. 
FigurelTUIshows the persistence function for the d = 3 FA 
model at different temperatures. Clearly, the dynamics 
is heterogeneous, and the spatial correlations of the local 
dynamics grow as T is decreased. The 'critical' nature 
of dynamic clusters is apparent: the pictures are remi- 
niscent of the spatial fluctuations of an order parameter 
close to a continuous phase transition, such as the mag- 
netization of an Ising model near criticality. In our case, 
the order parameter is a dynamic object, the persistence 
function, and the critical fluctuations are purely dynam- 
ical in origin |46j . 

We now quantify these observations. We can measure 
spatial correlations of the local dynamics via a spatial 
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FIG. 10: Spatial distribution of the local persistence at time 
t* such that P(t*) — 1/2 (i.e., 50% of sites, shown in black, 
have flipped by time t*). From top to bottom, T = 1.0, 0.2 
and 0.12. The appearance of dynamic critical fluctuations 
when T — > is evident. 



correlator of the persistence function, 

C(r,t,T) = j±- [(m)P i+r (t)) P 2 (t)] , (40) 

where the function f(t) = P(t) — P 2 (t) in the denomi- 
nator ensures the normalization C(r — 0,t.T) — 1. Al- 




FIG. 11: Time dependence of the dynamic susceptibility 14111 
at different temperatures. From left to right T = 1.0, 0.6, 0.4, 
0.3, 0.25, 0.2, 0.17, 0.15 and 0.13. The horizontal dotted line 
denotes the infinite time value, \{t — ► oo) = 1; the diagonal 
dotted line denotes the power-law fit x( r ,T) ~ r ' 46 . 



ternatively, one can take the Fourier transform of (I4U|) . 
giving the corresponding structure factor of the dynamic 
heterogeneity, 



S(q,t,T) 



Nf(t) 



^[(P fc (i)P/W)-P 2 W 



Jq(k-l) 



Finally, the zero wavevector limit of S(q, t) defines a dy- 
namic susceptibility, T) = S(q = 0,i, T), which can 
be rewritten as the normalized variance of the (unaver- 
aged) persistence function, pit) = N~ x Pi{t): 



x ( t ,T) = -^- ) [(p 2 (t))-(p(t)y< 



(41) 



Figure ^] shows the time dependence of the suscepti- 
bility H41fl for various temperatures. The behaviour of 
X is similar to that observed in atomistic simulations of 
supercooled liquids in general 0, and strong liquids in 
particular |6j. The susceptibility develops at low tem- 
perature a peak whose amplitude increases, and whose 
position shifts to larger times as T decreases. As ex- 
pected, the location of the peak scales with the relax- 
ation time t(T), indicating that dynamical trajectories 
are maximally heterogeneous when t ~ t(T). 

In Figure IT^l we show the correlator C(r,t,T) and the 
structure factor S(q 7 1, T) for different temperatures and 
fixed times t — t(T) where dynamic heterogeneity is 
maximal. These correlation functions clearly confirm the 
impression given by Fig. llUI that a dynamic length scale 
associated with spatial correlations of mobility develops 
and grows as T decreases. Note that at the lowest tem- 
peratures the structure factor does not reach a plateau 
at low q. This because the system size we use, although 
very large (N = 160 3 ), is not sufficiently so to allow us 
to probe the regime The necessary system sizes 

are simply too large to simulate on such long time scales. 

We can extract numerically the value of the dynamic 
length scale, £(T), at each temperature. To do so, we 
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(1 

FIG. 12: The spatial correlation function of dynamic hetero- 
geneity in real (top panel) and Fourier (bottom panel) spaces 
reveals the growth of a dynamic length scale as T — > 0. Tem- 
perature decreases from left to right (top panel) and from 
bottom to top (bottom panel). The dashed line in the bottom 
panel denotes the asymptotic behaviour S(q,r,T) ~ q~^ 2 ~ v \ 
when g£ S> 1, with n — —0.15. 

study in detail the shape of the correlation functions 
shown in Fig. 1121 As for standard critical phenomena, 
we find that the dynamic structure factor can be rescaled 
according to 

S(q,t,T)~ X (T,T)S(qO, (42) 

where the scaling function S(x) behaves as 

S(x ->■ 0) ~ const (43) 
S(x oo) ~ x 2 - 7 K (44) 

Both the susceptibility \ an d the dynamic length scale 
£ estimated at time t = r behave as power laws of the 
defect concentration, 

X~c-\ £~c-^. (45) 

These relations imply that the exponents 7 and v\_ should 
be numerically accessible by adjusting their values so that 
a plot of c 1 S versus qc~ VA - is independent of temperature. 
We show such a plot in Fig. ^3 and we find that the 
values 7 w 0.97 and v± rj 0.5 lead to a good collapse 




FIG. 13: Collapse of the dynamic structure factor of Fig. 1121 
using scaling laws 1451 and taking 7 = 0.97 and v± = 0.5. 
Dashed line is the scaling function S(x) = 1/(1 +X 2-77 ) with 
77 = -0.15. 



of the data. The exponent 7 can be independently and 
more directly estimated from Fig. 1111 by measuring the 
height of the maximum of the susceptibility for various 
concentrations. Fitting the result to a power law of c 
gives 7 ~ 0.96, in reasonable agreement with the first 
value. We find also that the scaling function S(x) is well- 
described by an empirical form S(x) = 1/(1 + x 2 ^ 71 ), 
consistent with Eq. (|44|) . Thus we can determine the 
value of the 'anomalous' exponent, 77; we find 77 w —0.15. 

As usual, it is difficult to estimate what constitutes the 
'best' collapse of the data, and so determine accurately 
the errors in the values of the exponents. Consequently, 
we are unable to determine v± with sufficient accuracy to 
conclude that it agrees — or disagrees — with the d = 3 DP 
value, v± p w 0.58. It is also difficult to compare 7 with 
its corresponding DP value, because this would require 
one to know the anomalous exponent 77 characterizing 
spatial correlations of the persistence function. From a 
field theory perspective this is a formidable task. How- 
ever, from our numerics we have that rj w —0.15, and so 
from scaling arguments we find 7 = (2 — r))v± ~ 1.075. 
This estimates lies however on the 'wrong' side of the 
classical value j cl — 1 as compared to the numerical value 
obtained above. 

We must conclude that numerical uncertainties are too 
large, and deviations from classical behaviour too small 
to make quantitative comparisons between DP and nu- 
merical exponents for spatial correlations. Plus, as we 
discussed above, our data may be subject at very low 
temperature to small, but unknown, finite size effects. 

We are nonetheless satisfied that the naive estimate 
v± = 1/d = 1/3 27] that one gets by estimating the mean 
distance between defects is invalidated by our numerical 
results. 
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VII. CONCLUSIONS 

We have derived a field theory for a kinetically con- 
strained model with isotropic facilitation, exemplified by 
the FA model. We have studied the field theory via RG, 
and the lattice-based FA model via numerical simula- 
tions. Our central results, briefly summarised in Ref. |2l| . 
are the following. 

The RG treatment suggests that the low-T dynamics 
is dominated by a non-classical, zero-temperature criti- 
cal point, which in turn implies that correlation times, 
dynamic correlation lengths and susceptibilities exhibit 
the scaling behaviour 

r~c~ A , £~c-" x , x~c" 7 , (46) 

with A = l+vii. The Arrhenius behaviour of the equilib- 
rium concentration of excitations, c « e _1 / T , gives rise to 
Arrhenius behaviour of the dynamics through (|46[1 . For 
dimensions d > 2, the critical point is that of DP, while 
for d = 1 it is that of CDP. The upper critical dimension 
is d c = 4, so that for dimensions d > 4 the exponents 
take classical values. For d = d c — 4 the exponents are 
classical, augmented with the usual logarithmic correc- 
tions 36]. For the time and space exponents we have 

A « 3, 2.3, 2.1, 2 (d= 1,2,3,>4), (47) 
v± w 1,0.73, 0.58, 1/2. (48) 

We have also performed large-scale numerical simula- 
tions of the FA model, which confirm many of the field- 
theoretic predictions. The relaxation times of the FA 
model, Figs. and \7\ follow the scaling laws given by 
(|46|l and l|47() in all dimensions simulated (d = 1 to 6). 
The existence of an upper critical dimension at d c = 4 
is evident (see Fig. EJ). The dynamics is increasingly 
heterogeneous and correlated in space as temperature is 
decreased, as can be seen, for example, in pictures of 
the local persistence (Tig. I10|l. The structure factor for 
this dynamic heterogeneity field in d = 3 exhibits scale- 
invariance (Fig. I13fl . 

More extensive simulations are required in order to 
clarify two further points. The spatial exponent obtained 
from the numerics is v± ~ 0.5, but we were unable to 
establish whether this number agrees precisely with the 
d = 3 DP value of i/± = 0.58. We also caution the reader 
that there may exist, even for d > 2, a crossover from 
early-time CDP behaviour to intermediate-time DP be- 
haviour, as is the case for some systems in, ostensibly, 
the DP universality class [2^ . 



Our work shows that standard theoretical methods, 
such as the renormalization group, can be used to ana- 
lyze coarse-grained models of glass-forming supercooled 
liquids |46l l47| . It supports the view that the dynam- 
ics of glass-formers is in many respects similar to that of 
standard critical phenomena, such as reaction-diffusion 
systems [35|,^3- We have found, numerically and analyt- 
ically, that the FA model and its associated field theory 
possess a zero-temperature critical point, in agreement 
with results obtained by other means [23. Rigorous re- 
sults confirm the existence of a T = critical point in 
other kinetically contrained systems, such as the East 
model |49|. and an analogous maximal-density critical 
point in the Kob- Andersen model [5(j. Extending the 
field theory treatment to models of fragile glass-forming 
liquids, such as the East model [34| and its generaliza- 
tions I23, constitutes an interesting challenge. 

Finally, our results provide some insight into the phys- 
ical meaning of fragility, in the Angell se nse 3 ■ First, we 
have shown here and elsewehere [l|| I23, |23| that strong 
systems show fluctuation-dominated heterogeneous dy- 
namics, in a similar manner to fragile systems. This con- 
tradicts the popular view that ' cooper ativity', 'fragility', 
and 'heterogeneity' are different facets of the same con- 
cept 00. 

In our view, the difference between strong and fragile 
liquids is in the strength of fluctuation effects. For ex- 
ample, the breakdown of the Stokes-Einstein relation ob- 
served in fragile liquids [Ell 152^ should also be observed 
in strong ones [53], but the effect will be less striking. 
However, since strong systems such as the FA model are 
characterized by a constant dynamic exponent, we expect 
that typical length scales at the glass transition are typ- 
ically larger in strong glass-formers than in fragile ones. 
Detailed studies of dynamic heterogeneity in atomistic 
models of strong liquids should be able to test these pre- 
dictions 0, while experimental investigations of strong 
glass-formers would also be very welcome. 
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